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Abstract 



Monte-Carlo simulations on a variety of 2d percolating systems at crit- 

icality suggest that the excess number of clusters in finite systems over the 

bulk value ric is a universal quantity, dependent upon the system shape but 

independent of the lattice and percolation type. Values of Uc are found to 

high accuracy, and for bond percolation confirm the theoretical predictions 

of Temperley and Lieb, and Baxter, Temperley and Ashley, which we have 

evaluated explicitly in terms of simple algebraic numbers. Predictions for the 

fluctuations are also verified for the first time. 
PACS numbers(s): 64.60Ak, 05.70.Jk 
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Typeset using REVT^ 



The standard percolation model [1] involves the random occupation of sites or bonds 
of a regular lattice. At a critical occupation probability Pc, the mean size of clusters of 
adjacently occupied sites becomes infinite, signalling the percolation transition. The number 
n{p) of clusters per site, however, remains finite, and attains a value ric — n(pc) with critical 
behavior \p — Pc\^''°' as p — > Pc, where 2 — a — 8/3. Critical exponents like a are universal, 
having the same value for all systems of a given dimensionality no matter what the lattice 
or percolation type, essentially because they derive from the large-scale, fractal properties 
of the percolation clusters. 

On the other hand, Uc, like Pc, is a non-universal quantity which varies from system 
to system. This is because it depends upon the microscopic, lattice-level behavior of the 
systems. For some systems, ric is known theoretically [2,3], while for others, it must be de- 
termined by Monte-Carlo (MC) means. Accurate values of ric — a seemingly basic property 
of a percolating system — do not seem to be available for many common lattices. Thus we 
embarked on a project to determine it for a variety of systems by MC simulation. 

To obtain accurate values, we employed a very extensive set of simulations. They were 
carried out on square (SQ) and triangular (TR) lattices, with both site (S) and bond (B) 
percolation, and the honeycomb (HC) lattice, for B only. Various orientations and system 
boundaries were used as shown in Fig. 1. Periodic boundary conditions were used. Systems 
ranged from 16 x 16 to 512 x 512 in size, and between 10^ to 10^ samples were generated for 
each, which required a total of several months of workstation computer time. All but one 
set of simulations were carried out at the critical point Pc of the respective lattice. 

Because simulations must necessarily be done on finite systems, an essential part of 
such a project is to determine the nature of the finite-size corrections. The use of periodic 
boundary conditions reduces those corrections to sub-surface terms. In characterizing the 
remaining corrections, we found a surprising result: for systems of a given shape, the finite- 
size corrections to ric appear to follow a universal behavior. For a finite system of S sites, 
we found that the average density of clusters (number per site) behaves as 
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n = nc + - + ... (1) 

where 6 is a function of the system shape only and thus appears to be a universal quantity 
For example, the results for a system with a square boundary of dimension L x L are shown 
in Fig. 2, where we plot n{S) vs. 1/S. For this boundary shape we considered three systems: 
S-SQ-I (site percolation on a square lattice with orientation I shown in Fig. 1), B-SQ-I, and 
B-SQ-II, where II is a square lattice rotated by 45°. The equality of the slopes of the 
three curves means that b has the same value 0.8835 for all these systems. Higher-order 
corrections to (1) are small and not visible above the statistical errors of our simulations. 

Considering rectangular systems of dimension L x L', we found that b increases mono- 
tonically with increasing aspect ratio r = L'/L as shown in Table I. Here, we considered 
only one type of system for each aspect ratio. Besides B-SQ, we also considered B-TR with 
the boundary III. With Ill's rectangular boundary, a normal torus is produced when the 
periodic boundary conditions are applied, and the effective aspect ratio of the system we 
used is Vs. Note that, because b{r) — b{l/r), the value 6(1) pa 0.8835 is a minimum for 
rectangular systems. 

Besides rectangular systems, we also considered a 60° rhomboid system with an effective 
hehcal or "twisted" periodic boundary conditions (systems IV and V) . This type of boundary 
results when a triangular lattice is represented in the computer by a square lattice with one 
set of diagonals drawn in. The periodic boundary conditions are applied to the lattice in 
the squared-off orientation, and an effective helical twist is introduced in the torus when 
the lattice is transformed back to its true configuration with the triangles equilateral. For 
such a system with sides L x L, we find b ~ 0.878, somewhat smaller than the minimum 
value above for a rectangular boundaries. (Likewise, 0.878 should be the minimum value of 
b for all 60° parallelograms.) Note that the effective height-to- width of this system is \/3/2, 
although it is not really comparable to a rectangular system of a given aspect ratio because 
of the twist in the boundary condition. 

While we do not have a proof that b is universal, we can make the following heuristic 
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arguments to support that hypothesis. Eq. (1) can be written 

N^ricS + b (2) 

where N — nS is the total number of clusters in the system. This equation shows that h 
is effectively the excess number of clusters over the expected bulk value, tt-c-S'. As such, it 
represents large clusters whose size is of the order of the size of the system. While the number 
and distribution of small, "microscopic" clusters is non-universal, the distribution of these 
larger clusters is. One can also relate h to the average number w of clusters wrapping around 
the toroidal system, which too is a universal quantity. Consider taking one 2L x 2L system, 
unwrapping the periodic boundary conditions, cutting it into four L x L systems, and then 
re-applying the periodic boundary conditions to each of these. The number of clusters N(L) 
should roughly satisfy N{2L) —w — 4{N{L) —w), since the density of non- wrapping clusters 
will be about the same, and this formula is consistent with (2) with w — b (here, -S" = L^) . A 
similar argument shows that b should increase as the aspect ratio r = L'/L increases, since 
there are then more wraparound clusters, consistent with what we have observed. 

The values of ric that we found from plots like Fig. 2 are listed in Table II. The results 
for the bond percolation systems are consistent with theoretical predictions that have been 
made for these systems. For B-SQ, where Pc = 1/2, Temperley and Lieb [2] showed 



B-SQ _ i 
~ 2 



z=i 



(3) 

IJ.=7r/3 



where Z = 2 cos and 

1 r°° , /7ra\ , /cosha — cos2u\ , 

/i = — / sech — In — ^ da , (4) 

4/^7-00 ^2//^ V cosha-l ^ ' 

and reported the numerical estimate O.O98O7. In fact, we have found that this rather com- 
plicated integral result can be evaluated explicitly [4]. First, integrate (4) by parts to obtain 

^ 4sinV /•-tan-i(tanhf^)cothf 

h — / '^—T da (5) 

TT Jo cosh a — cos 2a 

and then differentiate to find: 



dh 

dfj, 



SVS r°° tan ( tanh sinh 
TT io (2 cosh a + 1)2 



a 

—da 



/i=7r/3 

27 /■'^ a 

q; cosh — csch 3q; da . (6) 
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Substituting a = —2 logz in the first integral above, we can reduce it to a form that can be 
evaluated with the help of Mathematica [5], yielding 2^/3 — 3. Likewise, the second integral 
in (6) can also be evaluated, and leads to the simple final result 

^B-SQ ^ 3a/3 -5 ^ Q^gg . (7) 

Further details are given in [4]. Our numerical results in Table II are consistent with this 
prediction to all significant figures. Likewise, for B-TR, where pf'^^ = 2sin(7r/18) is the 
solution to — 3p + 1 = [6], Baxter, Temperley and Ashley [3] showed that 
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where 



i,^^-r 'j^^^t^^t^^^^dx , (9) 

2 J-oo X sinh nx cosh (f)x 
and estimated Uc ~ 0.1118 by numerical integration and also from the 16-term series avail- 
able at that time. Again, we have explicitly evaluated this expression [7,8], with the result 

^f'^'^^T-^ ^0-111844 275 (10) 

which is in agreement with our observations. 

Making use of duality, we also studied the B-HC system at the same time as B-TR. 
Our value of Uc for this lattice given in Table II agrees with prediction n^~^^ — nf"^^ -|- 
{Pc'^^Y ~ 0.153 733 341, which follows from the results of [16] written in terms of clusters 
of wetted sites per site on the TR lattice. Note that for all bond percolation systems, we 
count the number of clusters of all wetted sites, which includes both clusters of connected 
sites and "null" clusters of isolated sites with no occupied bonds attached to them. Having 
the above exact results for ric was very useful for finding b to high accuracy for these systems. 



Previous numerical results are also consistent with our findings. For B-SQ, Nakanishi 
and Stanley [9] found 0.98075 by simulation. Domb and Pearce [10] found Uc — 0.0173(3) 
by series analysis, where here Uc is the number of bond clusters per bond. To put (7) in this 
form, one must subtract the concentration of isolated sites (1/16) and then divide by two, 
the number of bonds per site. This yields 

24 a/3 _ 41 

nf ^'^(per bond) = ^ 0.017788 106, (11) 

and shows that the results of [10] are slightly low. Surprisingly, we have found little compar- 
ison with or discussion of the results of [2] and [3] in the percolation literature — a notable 
exception being [11], where in fact we first became aware of [2]. 
For S-SQ, where Pc = 0.592 746 0(5) [12], we find 

^s-SQ ^ Q Q27 5931 ± 0.000 0003 , (12) 

but could not find any published values to compare this with. A simple Pade analysis [13] 
of the 25-th order series for n{p) given in [14], with no attempt made to account for the 
branch-point singularity at Pc, yields Uc — 0.02754(4), the error representing the variation 
for [N,D] = [12,13], [13,12],.... 

For S-TR lattice, where Pc = 1/2; we find 

nf ^-R = 0.017 6255 ± 0.000 0005 . (13) 

This result is consistent with the previous MC value 0.017630(2) of Margolina et al. [15] 
but not the early low-order series value 0.0168(2) [10]. As in the S-SQ case, no theoretical 
prediction for this quantity exists. Note that the value of nf"^^ is quite close to nf'^^ 
of (11), which is reasonable because the matching site lattice to the B-SQ system is quite 
similar to the S-TR system, both having coordination number 6 and Pc — 1/2. 
We also considered S-SQ at p = 0.5, well below Pc- Our MC results give 

71^-^^(0.5) = 0.065 770 3 ± 0.000 000 2 . (14) 
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This quantity relates to a simple question of graph theory [17]: how many clusters of I's exist 
on a matrix filled randomly with O's and I's of equal concentration? Here finite-size effects 
were observed in lattices from 32 x 32 to 512 x 512 in size (i.e., 6 = 0), presumably because 
the dimension of a typical cluster was smaller than 32. In this case, a straightforward Pade 
analysis of Conway and Guttmann's series [14] yields a result of nearly the same accuracy: 
0.065 769 6(6). 

The fluctuations in the number of clusters from sample to sample, determined by 
((AiV)2) = (iV^) - {Nf = Su^ for large S, are given in Table II. Like ric, (t^ is a function 
of the lattice and percolation type and thus non-universal, but still independent of the 
boundary condition and shape. The value 0.164 found for B-SQ, for example, was found for 
square systems with both orientations, rectangular systems of aspect ratio 2 and 4, and also 
a square system with open boundaries. The error (one standard deviation) in the measured 
value of n is (j{SNruns)~^^'^ where N^uns is the number of samples; this is how we obtained the 
error bars in our numerical results above. Note that previously, only rough measurements 
of fluctuations were made [18]. 

For the B-SQ case, the theory of [2] yields a prediction for the fluctuations of a quantity 
H = Co + Sc: where Cq is the number of components (clusters plus isolated sites — the 
same as N above) and Sq is the number of independent cycles in the system: 



The numerical results given in [2] (p. 280) imply the value 0.1962 — 0.0377 = O.I585 for this 
quantity. After many intermediate steps, we have evaluated the resulting integral expressions 
exphcitly and flnd [4] 



At Pc, Sg — Cg: but because of the cross term {CgSg) in ((Aif)^), the above result cannot 
be used to obtain a prediction for = {{ACg)"^) /S. However, one can show that Sq is 
equal to the number of components on the dual lattice Cq, which can be measured at the 




(15) 



{{AHf) _ 8^-25 18 
S ~ 2 ^ V 



0.157781182 . 



(16) 



7 



same time as clusters on Cq- (It also follows that H is equal to the total number of hulls 
that can be drawn on the system, minus 1 if there is wraparound in both directions [19].) 
Indeed, for B-SQ percolation on rectangular systems, we find that ((Ai7)^) = 0.1578(1), in 
perfect agreement with (15). If Cq and Cq were uncorrelated, then {{AHY)/S would equal 
2(7^ = 0.328. This implies that realizations with more clusters on the regular lattice tend to 
have fewer clusters on the dual lattice, and vice versa, as one might expect. The fluctuations 
of H have finite-size corrections that behave as L~^, while the fluctuations of clusters 
alone was found to follow an unusual behavior fa L'^-^^. 

We have also measured {{AHY) / S for the B-TR case, where now Cq are components 
on the B-HC lattice, and find 0.2207(1). We have not related this result to theory [3]; 
although we can evaluate the next derivative of I2 [7], we have not found the additional 
terms analogous to the last two terms in (8). 

Thus, we have shown that the theory of [2] can be applied to cluster fluctuations if one 
includes the dual-lattice clusters in the measurements. We have found surprisingly simple 
expressions for the integrals of the theories of [2,3], which we have also verified numerically 
to high precision. And we have presented evidence that the excess number of clusters is a 
universal quantity, much like the crossing probability that has been of much interest recently 
[19-23]. Recently, Aharony and Stauffer [24] have presented theoretical arguments for the 
universality of b, and furthermore, Kleban [25] has derived explicit expressions for b using 
conformal invariance methods, which support this universality hypothesis. 

The authors thank P. Kleban, D. Stauffer and A. Aharony for their comments. This 
material is based in part upon work supported by the US National Science Foundation 
under Grant No. DMR-9520700. 
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FIGURES 

FIG. 1. Different lattices and boundary sfiapes used in the MC studies. Systems I, II and III 

have rectangular boundaries and lead to normal tori when the periodic boundaries arc applied. 
IV and V have a rhomboid boundary, but the periodic boundary conditions are applied when the 
lattice is in a squared-off form, leading to a helical or twisted boundary. 

FIG. 2. A plot of n = N/S vs. 1/S for the L x L S-SQ-I (filled square), B-SQ-I (open square) 
and B-SQ-II (open diamond) systems at Pc- The curves are offset vertically for clarity; the actual 
intercepts are nf'^^ for the upper curves, and nf'^^ for the lower two curves. This plot illustrates 
that h (the slope) is the same for all three systems. 
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TABLES 

TABLE I. Measured values of the universal finite-size correction constant b = N—ricS, grouped 
by universality class (system shape and boundary condition). Numbers in parentheses show the 
errors in the last digit (s). 



system 


boundary 


b 


S-SQT 


square 


0.8832(3) 


B-SQ-I 




8838f5) 


B-SQ-II 


square 


0.8835(8) 


B-TR-III 


1 X rectangle 


0.946(2) 


B-SQ-II 


1x2 rectangle 


0.991(2) 


B-SQ-II 


1x4 rectangle 


1.512(3) 


S-TR-IV 


1x1 rhombus 


0.8783(8) 


B-TR-IV 


1x1 rhombus 


0.878(1) 


B-HC-V 


1x1 rhombus 


0.877(1) 
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TABLE II. MC results for cluster concentration ric and fluctuations cr^ for the various systems 
studied, grouped by system type, along with known values of Pc [6,12]. These quantities are 
non-universal, but independent of the system boundary, as seen here for the B-SQ and B-TR 
cases. 
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S-SQ-I 
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B-SQ-I 


square 


0.5 


0.098 0763(8) 


0.1644 


B-SQ-II 


square 


0.5 


0.0980765(10) 


0.1644 


B-SQ-II 


1x2 rectangle 


0.5 


0.098076 


0.1642 


B-SQ-II 


1x4 rectangle 


0.5 


0.098 1 


0.163 


B-TR-III 


1 X ^/3 rect. 


0.347296 


0.111846(2) 


0.183 


B-TR-IV 


1x1 rhombus 


0.347296 


0.111843(2) 


0.1827 


B-HC-V 


1x1 rhombus 


0.652 705 


0.153 735(2) 


0.267 
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Fig. 2 - Ziff, Finch, and Adamciiik 
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